Astronomy & Astrophysics manuscript no. 13564 


©ESO 2010 


January 25, 2010 





Metal production in M33: space and time variations 

Laura Magrini 1 , Letizia Stanghellini 2 , Edvige Corbelli 1 , Daniele Galli 1 , Eva Villaver 3 

INAF-Osservatorio Astrofisico di Arcetri, Largo E. Fermi 5, 1-50125 Firenze, Italy 

e-mail: laura--edvige--galli@arcetri . astro, it 

National Optical Astronomy Observatories, Tucson, AZ 85719, USA 

e-mail: lstanghellini@noao.edu 

Universidad Autonoma de Madrid, Departamento de Ffsica Teorica C-XI, 28049 Madrid, Spain 
e-mail: eva.villaver@uam.es 



o 

(N 
G 

in 

(N 



O 

U 

Or 

6 

CO 



(N 
> 

(N 
(N 

On 
O 



X 



ABSTRACT 

Context. Nearby galaxies are ideal places to study metallicity gradients in detail and their time evolution. 

Aims. We analyse the spatial distribution of metals in M33 using a new sample and the literature data on H n regions, and constrain a 
model of galactic chemical evolution with H n region and planetary nebula (PN) abundances. 

Methods. We consider chemical abundances of a new sample of H n regions complemented with previous data sets. We compared 
H ii region and PN abundances obtained with a common set of observations taken at MMT. With an updated theoretical model, we 
followed the time evolution of the baryonic components and chemical abundances in the disk of M33, assuming that the galaxy is 
accreting gas from an external reservoir. 

Results. From the sample of H n regions, we find that i) the 2D metallicity distribution has an off-centre peak located in the southern 
arm; if) the oxygen abundance gradients in the northern and southern sectors, as well as in the nearest and farthest sides, are identical 
within the uncertainties, with slopes around -0.03-4 dex kpc -1 ; Hi) bright giant H n regions have a steeper abundance gradient than 
the other H n regions; f'v) H n regions and PNe have O/H gradients very close within the errors; v) our updated evolutionary model is 
able to reproduce the new observational constraints, as well as the metallicity gradient and its evolution. 

Conclusions. Supported by a uniform sample of nebular spectroscopic observations, we conclude that i) the metallicity distribution in 
M33 is very complex, showing a central depression in metallicity probably due to observational bias; ii) the metallicity gradient in the 
disk of M33 has a slope of -0.037± 0.009 dex kpc -1 in the whole radial range up to ~8 kpc, and -0.044± 0.009 dex kpc -1 excluding 
the central kpc; Hi) there is little evolution in the slope with time from the epoch of PN progenitor formation to the present. 

Key words. Galaxies: abundances, evolution - Galaxies, individual: M33 



1. Introduction 

Understanding the chemical evolution of nearby galaxies has 
become more and more interesting in the recent past, given 
the wealth of new, detailed data available for their stellar and 
gaseous components. M33 (NGC 598) is an ideal laboratory 
for testing chemical evolution models. As a nearby blue star- 
forming galaxy, with a large angular size (optical size 53'x83', 
Holmberg 1 1 958t . and an intermediate inclination (z'=53°), it is 
one of the galaxies observed with the greatest resolution and 
sensitivity. Distance estimates range from 730 kpc (Brunthaler 
et al. 120051) to 964 kpc (Bonanos et al. 2006). In thi s pape r 
we adopt an average value of 840 kpc (Freedman et al. 1 1 99 li t. 
Many spectroscopic studies are available, aiming at determin- 
ing the physical and chemical properties of its H n regions (cf., 
e.g., Smith [19751 Kwitt er & Alter [19811 Vfl chez et al. [19551 
W illner & Nelson-Patel [20021 Crocke tt et al. [20061 Ma grini et 
al. l2007al Rosolowsky & Simon 120081 Rubin et al. l2008l . These 
H n regions trace the present interstellar medium (ISM) and their 
metallicity is a measure of the star formation rate (SFR) inte- 
grated over the whole galaxy lifetime. As a result, the metallicity 
of H ii regions gives interesting constraints to galactic chemical 
evolution models. 

In the past, studies of the metallicity gradient of H n regions 
in M33 have shown a very steep profile. First quantitative spec- 
troscopic studies were carried out by Smith ( 119751) . Kwitter & 



Alter (11981b and Vflchez et al. G988 ). Their observations, lim- 
ited to the brightest and largest H n regions, led to a radial oxy- 
gen gradient about -0.1 dex kpc -1 . Considering the observations 
by the above researchers, Garnett et al. ( 1997) obtained an over- 
all O/H gradient of -0.1 1+0.02 dex kpc" 1 . 

With increasing sample sizes, flatter gradients have been 
determined, in particular by i) Willner & Nelson-Patel (2002) 
who derive neon abundances of 25 H n regions from infrared 
lines, obtaining a Ne/H gradient of -0.034+0.015 dex kpc -1 ; it) 
Crockett et al. (2006) who derive even shallower gradients (- 
0.016+ 0.017 dex kpc" 1 for Ne/H and -0.012+0.011 dex kpc" 1 
for O/H) from optical spectra of 6 H n regions; Hi) Magrini et 
al. (I2007al ) who obtain an O/H gradient of 14 H n regions, lo- 
cated in the radial range from ~2 to ~7.2 kpc with a slope of 
-0.054 ± 0.01 1 dex kpc" 1 ; iv) Rosolowsky & Simon (120081 who 
observe the largest sample of H n regions, 61, finding a slope of 
-0.027 ± 0.012 dex kpc -1 from -0.2 to ~6 kpc; v) Rubin et al. 
( 2008 ) who obtain Ne/H and S/H gradients from Spitzer infrared 
spectra of -0.058+0.014 dex kpc 1 and -0.052+0.021 dex kpc 1 , 
respectively. 

Rosolowsky & Simon (2008) attribute the large discrepan- 
cies in different authors' determinations to an intrinsic scatter of 
about 0. 1 dex around the average gradient, but the scatter is un- 
explained by the abundance uncertainties. This kind of scatter 
is also seen in the Milky Way gradient from H n regions (e.g. 
Afflerbach et al. 1 19971 ). and it may be caused by metallicity flue- 
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tuations in the ISM and by the spiral arms. Thus a limited num- 
ber of observations, coupled with a significant metallicity scatter 
at a given radius, may produce widely varying results. In the case 
of a shallow gradient this effect is even stronger; for example, for 
an abundance gradient of -0.02-0.03 dex kpc~' in a galaxy with 
a radius of 10 kpc and a scatter of 0. 1 dex, one would need more 
than 30 H n regions to obtain a good estimate of the slope of 
the gradient (Bresolin et al. 2009). Thus, only a large number 
of measurements can overcome the uncertainties engendered by 
the intrinsic variance and relatively shallow gradient of M33. 

Our sample is composed of 48 H n regions. We derived 
the physical and chemical properties for 19 H n regions which 
have not been observed previously, and for other 14 H n regions 
whose chemical abundances have already been published. For 
the remaining H n regions, the faintness of their spectra did not 
allow any reliable abundance determination, since their electron 
temperature (T e ) could not be derived. We complemented these 
observations by a sample of 102 planetary nebulae (PNe), al- 
ready presented by Magrini et al. (2009a, hereafter M09), ob- 
served during the same run. The main advantage of observing 
a combined sample of H n regions and PNe is being able to use 
not only the same observational set-up, but also the same data re- 
duction and analysis techniques, and to use identical abundance 
determination methods. 

Although our sample of H n regions does not add much to the 
literature, the presence of several objects in common with pre- 
vious studies allows us to check the consistency of different sets 
of chemical abundance results. By including at the same time 
two stellar populations of different ages but with similar spec- 
troscopic characteristics, our observations allowed us to study 
for the first time the global metallicity, its 2-dimensional (2D) 
distribution and its radial gradient, at two different epochs in the 
galaxy's lifetime avoiding biases introduced by different metal- 
licity analysis. The aim of the present study is to settle the ques- 
tions of the value of the metallicity gradient in M33 and its time 
evolution. In this framework, the new observations of H n re- 
gions and PNe in M33 complemented with the previous data 
represent the largest metallicity database available for an exter- 
nal galaxy. 

In addition to metallicity data, recent results, such as the de- 
tection of inside-out growth in the disk of M33 (Williams et 
al. 12009). the detailed analysis of the star formation both in the 
whole disk (Verley et al. 2009) and in several giant H n regions 
(Relano & Kennicutt 2009), stimulated us to revise the already 
existing chemical evolution model (M07b) and the star forma- 
tion process in M33. Particular attention was put on the obser- 
vational constraints that our previous model failed to reproduce, 
such as the radial profile of the molecular gas and the relation- 
ship between the SFR and the molecular gas. 

The paper is organized as follows. In Sect. [2] we describe our 
sample of H n regions observed with MMT These data, together 
with a large literature dataset, allowed us to compute the metal- 
licity gradient of H n regions again. In Sect. [3] we present the 2D 
distribution of the metallicity and the radial gradient of different 
types of H n regions. In Sect. 4 we discuss the off-centre metal- 
licity peak and its origin. In Sect. [5] the data are compared with 
the prediction of chemical evolution model of M33. Finally, our 
conclusions and a summary are given in Sect. [6] 

2. The H ii region dataset 

Hot O-B stars ionize their surrounding medium, producing the 
characteristic emission-line spectra of H n regions. The H n re- 
gions of M33 studied in the literature span a wide range of lu- 



minosities. Their intrinsic brightness led giant H n regions to 
be preferred in the earlier studies (e.g. Kwitter & Aller 1981 
Vilchez et al. 1988) when only relatively small telescopes were 
available. Smaller and fainter H n regions have instead been the 
subject of later spectroscopic investigations (e.g., Magrini et al. 
2007a Rosolowsky & Simon 2008). The latest abundance deter- 
minations have been restricted to the H n regions with available 
electron temperature measurements. Several emission-line diag- 
nostics of nebular T e are indeed present in the optical spectrum 
of an H ii region, namely: [O m] 4363 A, [N n] 5755 A, [S in] 
6312 A, [O n] 7320-7330 A. Determining T e is the only way to 
derive the ionic and total chemical abundances safely and accu- 
rately. An assumed T e could produce error of a factor of 2 or 
more in the final chemical abundances (cf., e.g., Osterbrock & 
Ferland 1989). This is why in the following analysis we include 
only those H n regions whose T e is directly measured. 

2.1. The MMT observations: data reduction and analysis 

In November 2007, we obtained spectra of 48 H n regions (and 
102 PNe) in M33 usin g the M MT Hectospec fiber-fed spectro- 
graph (Fabricant et al. 2005). The spectrograph was equipped 
with an atmospheric dispersion corrector and it was used with a 
single setup: 270 mirT 1 grating at a dispersion of 1.2 A pixel 1 . 
The resulting total spectral coverage ranged from approximately 
3600 A to 9100 A, thus including the basic emission-lines neces- 
sary for determining their physical and chemical properties. The 
instrument deploys 300 fibers over a 1 -degree diameter field of 
view, and the fiber diameter is ~ 1.5" (6 pc adopting a distance 
of 840 kpc to M 33). 

Some of the H ii regions in our sample already have pub- 
lished spectra in the literature so we use them as control sample, 
while several are new. In Table Q] we list the H n regions from 
the total observed sample for which we can derive the physi- 
cal and chemical properties. In Table Q] we list the H ii regions 
from the total observed sample for which we can derive the 
physical and chemical properties. The identification names are 
from: BCLMP- Boulesteix et al . (fT974l ; CPSDP- Courtes et al. 
(fT987l>; G DK99-Gordon et al. (fT996l ; EPR2 003-Engargiola et 
al. (I20U31 ; MJ98-Massey & Johnson JT9981 . The H n regions 
not identified in previous works are labled with LGC-HII-n as in 
Magrini et al. (2007a), standing for H n regions discovered by 
the Local Group Census project (cf. Corradi & Magrini 2006). 
The coordinates J2000.0 of the position of the fibers projected 
on the sky are shown in the third and forth columns. They do not 
correspond exactly to the centre of the emission line objects, but 
generally to the maximum [O in] emissivity. 

The spectra were reduced using the Hectospec package. The 
relative flux calibration was done observing the standard star 
Hiltm600 (Massey et al. 119881 ) during the nights of October 15 
and November 27. The emission-line fluxes were measured with 
the package SPLOT of IRAlfl Errors in the fluxes were calcu- 
lated taking the statistical error in the measurement of the fluxes 
into account, as well as systematic errors of the flux calibrations, 
background determination, and sky subtraction. The observed 
line fluxes were corrected for the effect of the interstellar ex- 
tinction using the extinction law of Mathis ( 1990) with Ry-3.1. 
We derived c(H/3), the logarithmic nebular extinction, by using 
the weighted average of the observed-to-theoretical Balmer ra- 

1 IRAF is distributed by the National Optical Astronomy 
Observatory, which is operated by the Association of Universities 
for Research in Astronomy (AURA) under cooperative agreement with 
the National Science Foundation 
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Table 1. MMT observations of H n re 
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dance determination. 








New sample 


ID 


tt"a 

RA 


Dec 




(1) 


(2) 


(3) 


(4) 




1 


BCLMP 275 A 




30:36:07.90 




2 


GDK99 3 


1:32:31.7 


30:35:27.39 




3 


LGCHII14 


1:32:33.3 


30:32:01.90 




4 


CPSDP 26 


1:32:33.7 


30:27:06.60 




5 


LGCHII15 


l:32:4U.fi 


30:24:24.99 




6 


EPR2003 87 


1 . JZ.4Z.4 


30:22:25.59 




7 


LGCHII16 


1:32:43.5 


30:35:17.29 




8 


LGCHII17 


1:33: 11.3 


30:39:03.39 


1 


9 


BCLMP 694 


1:33:52.1 


30:47:15.59 


10 


BCLMP 759 


1:33:56.8 


30:22:16.50 


E 

o 


1 1 


MJ98 WR 112 


1:33:57.3 


30:35:11.09 


12 


LGCHII18 


1:33:58.9 


30:55:31.30 


QC 


13 


BCLMP 282 


1:32:39.1 


30:40:42.10 


U 

aj 


14 


BCLMP 264 


1:32:40.2 


30:22:34.70 




15 


BCLMP 238 


1:32:44.5 


30:34:54.30 


y, 


16 


BCLMP 239 


1:32:51.8 


30:33:05.20 


Flu 


17 


BCLMP 261 


1:32:54.1 


30:23:18.70 


18 


CPSDP 123 


1:33:20.4 


30:32:49.20 




19 


CPSDP 43A 


1:33:23.9 


30:26:15.00 





Control sample 



20 


LGCHII2 


1:32:43.0 


30:19:31.19 


21 


LGCHII3 


1:32:45.9 


30:41:35.50 


22 


BCLMP289 


1:32:58.5 


30:44:28.60 


23 


BCLMP218 


1:33:00.3 


30:30:47.30 


24 


MAI 


1:33:03.4 


30:11:18.70 


25 


BCLMP290 


1:33:11.4 


30:45:15.09 


26 


IC132 


1:33:15.8 


30:56:45.00 


27 


BCLMP45 


1:33:29.0 


30:40:24.79 


28 


BCLMP670 


1:34:03.3 


30:53:09.29 


29 


MA2 


1:34:15.5 


30:37:11.00 


30 


BCLMP691 


1:34:16.6 


30:51:53.99 


31 


IC131 


1:33:15.0 


30:45:09.00 


32 


IC133 


1:33:15.9 


30:53:01.00 


33 


BCLMP745 


1:34:37.6 


30:34:55.00 



BCLMP 694 R=2 kpc 



-JO acc 



MAI R= 7.4 kpc 



m SA o 



wavelength (A) 

Fig. 1. Two example spectra of H n regions located at different 
galactocentric distance: BCLMP 694 at ~ 2 kpc and MAI at 
-7.4 kpc. 



tios of Ha, Hy, and H<5 to U/3 (Osterbrock & Ferland [1989l . The 
detailed description of the data reduction and the plasma and 
chemical analysis can be found in Magrini et al. (2009a, here- 
after M09). Spectra of two H n regions, one close to the galactic 
centre and one in the outer part of the disk are shown in FigJTJ 

Table [2] gives the results of our emission-line measurements 
and extinction corrections for 33 H n regions whose spectra were 
suitable for determining physical and chemical properties. The 
columns of Table [2] indicate: (1) H n region name; (2) nebular 
extinction coefficient c(H/3) with its error; (3) emitting ion; (4) 
rest-frame wavelength in A; (5) measured line fluxes; (6) abso- 
lute errors on the measured line fluxes; (7) extinction corrected 
line fluxes. Both (5) and 1^ (7) are expressed on a scale where 
Hy8=100. Table|2]is published in its entirety in the electronic edi- 
tion of A&A. A portion is shown here for guidance regarding 
its form and content. The analysed H n regions represent about 
2/3 of our sample. The remaining 1/3 H n regions have noisy 
spectra and are distributed at all galactocentric radii. We used 
the extinction-corrected intensities to obtain the electron densi- 
ties and temperatures. Electronic density was derived from the 
intensities of the sulphur-line doublet [S n] 6716,6731 A. We 
used the intensities of several emission-line ratios, when avail- 
able, to derive low and medium-excitation temperatures (see also 
Osterbrock & Ferland [T989l §5.2): [N ii]A5755/(/16548 + ,16584) 
and [O n]/13727/(i7320 + ,17330) for low-excitation T e , while 
[O m],14363/(,15007 + ,14959) and [S m] ^63 12/(^9069 + ,19532) 
for medium-excitation T e . We performed plasma diagnostics and 
ionic abundance calculation by using the 5-level atom model in- 
cluded in the nebular analysis package in IRAF/STSDAS (Shaw 
& Dufour 1994). The elemental abundances are then determined 
by applying the ionization correction factors (ICFs) following 
the prescriptions by Kingsburgh & Barlow d 19941 ) for the case 



where only optical lines are available. In the abundance analysis 
we adopted T e [N II] and/or T e fO II] for computing the N + , + , 
S + abundances, while T e [0 III] and/or T e [S III] for 2+ , S 2+ , 
Ar 2+ , He + , and He 2+ . We calculated the abundances of He i and 
He H using the equations of Benjamin et al. ( 1999) in two den- 
sity regimes, i.e. n e >1000 cirT 3 and <1000 crrT 3 . The Clegg's 
collisional populations were taken into account (Clegg 1987J. In 
Table [3] we present the electron densities and temperature, and 
the ionic and total chemical abundances of our H n region sam- 
ple, which only includes H n regions with at least one measured 
value of T e . The columns of Table [3] present: (1) identification 
name; (2) label of each plasma diagnostic and abundances avail- 
able; (3) relative values obtained from our analysis. Table 3 is 
published entirety in the electronic edition. 

We derived the temperature and density uncertainties using 
the error propagation of the absolute errors on the line fluxes. 
The errors on the ionic and total abundances were computed 
taking the uncertainties in the observed fluxes, in the electron 
temperatures and densities, and in c(Hy6)into account. In Table |4] 
a summary of the total abundances He/H, O/H, N/H, Ne/H, S/H. 
and Ar/H and their errors are presented. The He abundance is 
shown by number with its absolute error, while the metal abun- 
dances are expressed in the form of log(X/H)+12 with errors 
expressed in dex. The last row indicates the average abundances 
computed by number. 

2.2. The PN data-set 

The PN population of M33 was studied by Magrini et al. (f2009a) 
using multi-fiber spectroscopy with Hectospec at the MMT with 
the same spectroscopic setup and during the same nights as the 
H ii region observations presented here. Spectra of 102 PNe were 
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Table 2. Observed and de-reddened fluxes. 



ID c(H/3) 


Ion 


A (A) 


Fj 


A(Fj) 


h 


(1) (2) 


(3) 


(4) 


(5) 


(6) 


(7) 


BCLMP275A 0.493±0.009 


[Oil] 


3727 


186.4 


2.2 


256.7 




[NeIII]/HI 


3968 


3.2 


0.7 


4.1 




HI 


4100 


13.0 


1.1 


16.2 




HI 


4340 


34.4 


1 . 1 


39.8 




[OIII] 


4363 


1.1 


0.8 


1 .2 




HI 


4861 


100.0 


1.6 


100.1 




[OIII] 


4959 


57.3 


1 .2 


55.8 




[OIII] 


5007 


169.9 


2.0 


163.3 




Hel 


5876 


13.5 


0.9 


10.7 




[SHI] 


63 1 2 


1 .7 


0.8 


1 .3 




[Nil] 


6548 


9.1 


1 .0 


6.5 




HI 


6563 


409.6 


2.3 


290.2 




[Nil] 


6584 


28.9 


1 . 1 


20.4 




Hel 


6678 


3.9 


0.8 


2.7 




[SH] 


6717 


32.7 


1 . 1 


22.6 




[SII] 


6731 


23. 1 


1 .2 


15.9 




Hel 


7065 


2.6 


0.8 


1 .7 




[Arlll] 


7135 


9.8 


0.8 


6.4 




[srn] 


9069 


21.4 


1.0 


10.2 


GDK99 3 0.448±0.007 


[Oil] 


3727 


83.3 


1.0 


1 1 1.4 




HI 


3835 


6.8 


0.6 


8.8 




[Nelll] 


3869 


25.9 


0.6 


33.4 




Hel 


3889 


9.8 


0.6 


12.6 




[NeIII]/HI 


3968 


18.0 


0.7 


22.7 




HI 


4100 


17.5 


0.7 


21.2 




HI 


4340 


40.0 


0.9 


45.8 




[OIII] 


4363 


3.3 


0.5 


3.8 




Hel 


4471 


3.5 


0.5 


3.8 




Hell 


4686 


0.5 


0.5 


0.6 




HI 


4861 


100.0 


1 .3 


100.1 




[OIII] 


4959 


196.2 


1.5 


191.5 




[OIII] 


5007 


585.9 


2.8 


565.0 




Hel 


5876 


14.0 


0.7 


11.3 




[Nil] 


6548 


6.1 


0.4 


4.5 




HI 


6563 


396.4 


1.9 


289.9 




[Nil] 


6584 


19.2 


0.6 


14.0 




Hel 


6678 


4.3 


0.4 


3.1 




[SII] 


6717 


31.0 


0.7 


22.2 




[SII] 


6731 


21.6 


0.7 


15.4 




Hel 


7065 


3.4 


0.5 


2.3 




[Arlll] 


7135 


14.8 


0.5 


9.9 



analysed and plasma diagnostics and chemical abundances ob- 
tained for 93 PNe where the necessary diagnostic lines were 
measured. The data reduction and the plasma diagnostics fol- 
lowed exactly the same procedure as described in the present 
paper, thus ensuring that no biases are introduced for the differ- 
ent analysis of the spectra. About 20% of the studied PNe have 
young progenitors, the so-called Type I PNe. The rest of the PNe 
in the sample are the progenies of an old disk stellar population, 
with main sequence masses M<3M Q and ages t>0.3 Gyr. A tight 
relation between the O/H and Ne/H abundances was found, ex- 
cluding that both elements have been altered by PN progenitors 
and supporting the validity of oxygen as a good tracer of the ISM 
composition at the epoch of the progenitors' birth. 



3. The metallicity distribution in M33 

The large amount of chemical abundance data from H n regions 
in M33 allow us to analyse the spatially-resolved distribution of 
metals in the ISM. In this section, we present the radial distri- 
bution and the map of O/H, using the new data presented in this 
paper and all previous oxygen determinations for which T e has 
been measured. 



Table 3. Plasma diagnostics and abundances. 



ID 



(1) 


(2) 


(3) 


BCLMP725A 








T e [0 III] 


10600 




T [S III] 


15800 




Hel/H 


0.076 




[OII]/H 


7.740e-05 




[01II]/H 


4.737e-05 




ICF(O) 


1.000 




O/H 


1.248e-04 




[NII]/H 


3.287e-06 




ICF(N) 


1.612 




N/H 


5.299e-06 




[ArUrj/H 


4.910e-07 




ICF(Ar) 


1.87 




Ar/H 


9.182e-07 




[SII]/H 


7.486e-07 




[SIII]/H 


1.100e-06 




ICF(S) 


1.019 




S/H 


1.884e-06 


GDK99 3 








T e [OIII] 


10200 




Hel/H 


0.081 




[OII]/H 


3.880e-05 




[01II]/H 


1.860e-04 




ICF(O) 


1.006 




O/H 


2.261e-04 




[NII]/H 


2.536e-06 




ICF(N) 


5.826 




N/H 


1.475e-05 




[NeIII]/H 


3.150e-05 




ICF(Ne) 


1.215 




Ne/H 


3.829e-05 




[ArUTJ/H 


8.560e-07 




ICF(Ar) 


1.87 




Ar/H 


1.601e-06 




[SII]/H 


8.157e-07 




ICF(S) 


1.323 




S/H 


8.255e-06 



3.1. The metallicity gradient of H n regions 

Our cumulative sample includes: i) H n regions by Magrini et al. 
(2007a ), which includes abundances from their own sample and 
previous abundance determinations (all with T e and with abun- 
dances recomputed uniformly); ii) the sample by Rosolowsky & 
Simon (2008 ); Hi) the present sample (see Table [3). In Fig. [2] 
we show the oxygen abundance as a function of galactocentric 
distance for the cumulative sample of H n regions. In this figure 
each point corresponds to a single region; i.e., we do not plot 
multiple measurements for the same region but only the value 
with the lowest error. Note the large dispersion in the radial re- 
gion between 1 and 2 kpc from the centre caused by several 
high- and low-metallicity regions, located in the southern arm 
(see Sect. 3.3), which might be related to the presence of a bar 
(e. g., Corbelli & Walterbos 2007 ). We applied the routine ^zfexy 
in Numerical Recipes (Press et al. 1992) to fit the relation be- 
tween the oxygen abundances and the galactocentric distances, 
taking their errors into account and minimizing^ 2 . Typical er- 
rors on the de-projected galactocentric distances associated with 
the uncertainty on the inclination were less than 0. 1 kpc (Magrini 
et al. [2007a| . The fit to the complete sample gives a gradient of 

12 + log(0/H) = -0.037(+0.009) R GC + 8.465(+0.038) (1) 

where Rgc is the de-projected galactocentric distance in kpc, 
computed by assuming an inclination of 53° and a position an- 
gle of 22°. A weighted linear least-square fit to the MMT sample 
only gives a gradient of 

12 + log(0/H) = -0.044(+0.017) R GC + 8.447(+0.084), (2) 

which is consistent within the errors with the gradient from the 
larger sample. In the rest of the paper, we use the larger and more 
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Table 4. The chemical abundances of our MMT sample. 



# 


He/H 


O/H 


N/H 


Ne/H 


Ar/H 


S/H 


(1) 


(2) 


(3) 


(4) 


(5) 


(6) 


(7) 


1 


0.076±0.001 


8.10±0.07 


6.72±0.15 


- 


5.96±0.34 


6.27±0.20 


2 


0.081 ±0.001 


8.35±0.05 


7.17±0.11 


7.58±0.015 


6.20±0.09 


6.92±0.21 


3 


0.058±0.003 


7.97±0.08 


6.88±0.12 


- 


6.02±0.11 


6.62±0.15 


4 


0.088±0.010 


8.36±0.05 


7.24±0.19 


- 


5.82±0.18 


6.46±0.19 


5 


0.085±0.005 


8. 45 ±0.02 


7.35±0.14 


- 


- 


6.59±0.12 


6 


0.114±0.003 


8.52±0.08 


6.86±0.15 


7.60±0.20 


6.27±0.12 


6.50±0.20 


7 


0.085±0.005 


8.30±0.12 


7.09±0.20 


7.00±0.18 


6.07±0.16 


6.61±0.15 


8 


0.113±0.005 


8.21±0.07 


7.30±0.24 


7.69±0.42 


6.20±0.20 


6.49±0.15 


9 


0.092±0.008 


8.33±0.10 


7.21±0.10 


7.53±0.20 


6.22±0.30 


6.79±0.10 


10 


0.077±0.005 


8. 11 ±0.08 


6.88±0.29 


- 


6.16±0.22 


6.66±0.15 


1 1 


0.093±0.003 


8.30±0.07 


7.37±0.02 


- 


6.30±0.3O 


7.09±0.11 


12 


0.091 ±0.008 


8. 15 ±0.07 


7.94±0.05 


7.20±0.10 


6.24±0.10 


6.51 ±0.20 


13 


0.080±0.005 


8.32±0.05 


6.92±0.12 


- 


6.19±0.12 


6.36±0.18 


14 


0.076±0.003 


7.87±0.08 


6.93±0.15 


7.11±0.20 


5.89±0.20 


6.41±0.18 


15 


0.089±0.006 


7.92±0.08 


6.75±0.15 


- 


5.92±0.18 


6.57±0.15 


16 


0.086± 0.005 


8.27±0.07 


7.20±0.14 


- 


- 


6.52±0.20 


17 


0.115±0.008 


8.01±0.08 


6.62±0.16 


7.16±0.18 


5.93±0.14 


6.40±0.18 


IS 


0.078±0.005 


8.35±0.05 


7.35±0.12 


- 


6.35±0.15 


6.94±0.20 


19 


0.059±0.003 


8.40±0.05 


8.05±0.10 


8.03±0.18 


5.16±0.18 


7.02±0.17 


20 


0.083±0.001 


8.08±0.05 


6.91±0.15 


6.87±0.18 


6.18±0.25 


6.64±0.15 


21 


0.086±0.001 


8.42±0.06 


7.56±0.13 


7.31±0.15 


6.01±0.27 


6.63±0.10 


22 


0.072±0.005 


8.35±0.12 


7.34±0.16 


7.75±0.38 


5.67±0.50 


6.89±0.11 


23 


0.088±0.001 


8.17±0.12 


6.97±0.18 


- 


6.15±0.14 


6.73±0.23 


24 


0.080±0.008 


8.28±0.15 


7.10±0.20 


7.61±0.28 


6.14±0.32 


6.71±0.40 


25 


0.096±0.005 


8.38±0.13 


7.37±0.15 


7.36±0.20 


5.84±0.10 


6.57±0.15 


26 


A CU< T AA"2 


/.ya±o.o5 


o.yfi±o. 15 


H 10 i A 11 
/.Z5±U. 12 


5.86±0. 15 


6.36±0. 13 


27 


0.095±0.001 


8.48±0.08 


7.62±0.12 


7.73±0.15 


6.33±0.25 


6.87±0.15 


28 


0.088±0.002 


8.30±0.07 


7.10±0.18 


7.42±0.20 


6.23±0.30 


6.63±0.20 


29 


0.091 ±0.005 


8.31±0.10 


7.19±0.15 


7.46±0.25 


6.42±0.24 


6.92±0.21 


30 


0.096±0.003 


8.42±0.06 


7.20±0.12 


7.85±0.21 


6.37±0.32 


6.75±0.23 


31 


0.097±0.005 


8.47±0.08 


7.26±0.15 


7.76±0.20 


6.29±0.25 


7.06±0.15 


32 


0.079±0.005 


8.27±0.08 


7.21±0.17 


7.58± 0.21 


5.49± 0.30 


6.83±0.15 


33 


0.067±0.008 


7.93±0.10 


7.10±0.20 




6.09±0.21 


6.51±0.20 




0.085±0.01 1 


8.27!^ 


n ? i +0.30 
'■ J1 -0.35 


7 56 +0JS 


6 13+ 014 


6 71 +019 

°'' 1 _0.34 



complete sample when discussing the metallicity gradient and its 
possible time variation, but excluding the first kpc region where 
only a few low metallicity regions were analysed. We discuss 
the possible reasons for the lower value of the central metallicity 
later in this section and in Sect. [4] The O/H gradient of the whole 
sample of H n regions sample, excluding the central 1 kpc, is 

12 + log(0/H) = -0.044(+0.009) R GC + 8.498(+0.041). (3) 

We also checked the metallicity gradient by tracing it in dif- 
ferent areas of the galaxy, namely in the northwest and in the 
southeast halves, separately, and in the nearest and in the far- 
thest sides. The results are shown in Fig. [3] The northern and 
southern gradients, as well as those relative to the nearest and 
farthest sides, are identical within the uncertainties, with slopes 
around -0.03-4 dex kpc 1 . The only difference found between 
the metallicity gradients obtained for sections of the galaxy is 
the presence of a high metallicity peak in the southern arm. In 
conclusion, the present H n region sample (literature + present- 
work, 103 objects), including only nebulae with measured T e , 
reinforce the recent results on the slope of the M33 O/H gradi- 
ent, with a global slope up to around 8 kpc of -0.03 dex kpc -1 , 
and of about -0.04 dex kpc -1 excluding the central 1 kpc. The 
very central regions remain somewhat undersampled (6 objects 
within 1 kpc from the centre, and 9 within 1.5 kpc) and in dis- 
agreement with other results. A comparison with the metallicity 
gradient derived from young stars, which are representative of 
the same epoch in the lifetime of the galaxy as H n regions, and 
with the infrared spectroscopy of H n regions, is necessary. 

Stellar abundances were obtained by Herrero et al. ( 1994) for 
AB-supergiants, McCarthy et al. d 19951 1 and Venn et al. d 19981 
for A-type supergiant stars, and Monteverde et al. ( 119971 [2000) 
and Urbaneja et al. ( 120051 1 for B-type supergiant stars. The 



largest sample of Urbaneja et al. (2005]) gave a O/H gradient 
of -0.06+0.02 dex kpc -1 . Recently, U et al. ( [20091 has presented 
spectroscopic observations of a set of A and B supergiants. They 
determined stellar metallicities and derived the metallicity gra- 
dient in the disk of M33, finding solar metallicity at the centre 
and 0.3 solar in the outskirts at a distance of 8 kpc. Their average 
metallicity gradient is -0.07+0.01 dex kpc -1 . At a given radius, 
H n regions have abundances slightly below the stellar results, 
and this is probably due to the depletion of oxygen in H n re- 
gions on dust grains (e.g., Bresolin et al. [2009). The slopes of 
the two gradients agree if the comparison is done between about 
1 ant 8 kpc. The cause of the difference between the supergiant 
and H n region gradient is the metallicity value in the central re- 
gions. In fact, the H n regions located within 1 kpc from centre 
have metallicity below solar, whereas the supergiants are metal 
rich, ranging from solar values to above solar. The origin of this 
discrepancy is not the temperature gradients within the nebulae 
(Stasiriska 2005) because they become important at higher metal 
abundances. 

Recent observations of 25 H n regions by Rubin et al. (120081 ) 
with Spitzer have allowed a measurement of the Ne/H and S/H 
gradient across the disk of M33 showing no decrease in chemical 
abundances in the central regions. Infrared Ne and S emission 
lines do not have a strong dependence on T e , and consequently 
their abundances can be determined even without a temperature 
measurement. 

One way to explain the low metallicity in the central 1.0x1.0 
kpc 2 area is related to the criterion used for H n regions. Usually, 
chemical abundances derived from optical spectroscopy rely on 
direct measurement of the electron temperature, given by the 
[O m] 4363A emission line. This emission line is inversely 
proportional to the oxygen abundance and barely detectable for 
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Fig. 2. The O/H radial gradient for the cumulative H n region 
sample: filled circles are the MMT observations (new and con- 
trol samples), empty circles are the literature abundances. The 
continuous line is the weighted linear least-square fit of Equation 
3, i.e., with a radial range from 1 to 8 kpc from the M33 centre. 
The dashed vertical line indicates the regions located at its left- 
side and excluded from the fit. 



0/H>8.6 from an average luminosity nebula (e.g., Nagao et al. 
2006 ). The request for [O m] 4363 A detection might determine 
a bias towards lower metallicity, with the exclusion of the highest 
metallicity nebulae. This could explain the differences between 
the optical spectroscopy results and both the stellar abundance 
determinations, and the H n region infrared spectroscopy. 

3.2. The abundance gradients of the other chemical 
elements 

Our analysis allowed us to measure other chemical elements in 
addition to oxygen, as He/H, N/H, Ne/H, S/H. and Ar/H. The 
more reliable measurement is that of oxygen for the reasons il- 
lustrated in the Appendix, and we use it to follow the chemical 
evolution of M33. Nevertheless the other chemical elements are 
measured in enough H n regions to compute their radial gradi- 
ents. 

In Table [5] we show the slopes and central abundance val- 
ues of the radial gradients of N/H, Ne/H, S/H. and Ar/H of our 
sample of H n regions. We did not calculate the radial gradi- 
ent of He/H because we measured only the ionized fraction of 
He in H n regions, which is only a small part of the total he- 
lium abundance. All gradients have a negative slope, consistent, 
within the errors, with the slope found for O/H, while N/H is 
a bit steeper, as already noticed, e.g., by Magrini et al.(fJ007a). 
Its different behaviour with respect to the <r-elements, as oxygen, 
neon, sulphur, and argon, comes from the different places of pro- 
duction, a-elements are indeed produced by SNe II, which are 
the final phase of the evolution of massive stars, while nitrogen 
is one of the final products of the evolution of long-lived low- 
and intermediate-mass stars. This is discussed in detail in Sec|5] 
Finally, there is a very good agreement of the S/H and Ne/H 
gradients with those derived from the infrared spectra of H n re- 
gions by Rubin et al. (120081 1. for which they found a gradient of 
-0.058+0.014 dex kpc" 1 for Ne/H and -0.052+0.021 dex kpc" 1 
for S/H. 




013345678 
R (kpc) 



Fig. 3. The O/H radial gradient obtained in different regions of 
M33. Top: nearest side (filled squares) and farthest side (empty 
squares) gradient. Bottom: North (filled squares) and South 
(empty squares) gradient. In each panel the weighted linear least- 
square fits of the two regions are shown with two lines (continu- 
ous and dotted). 

Table 5. The radial gradients of N/H, Ne/H, S/H. and Ar/H 



12 + loa(X/H) 


slope 


central value 


(1) 


(2) 


(3) 


N/H 


-0.08+0.03 


7.53+0.15 


Ne/H 


-0.05+0.04 


7.71+0.21 


S/H 


-0.06+0.02 


6.41+0.11 


Ar/H 


-0.07+0.03 


6.98+0.13 



3.3. The population-dependent metallicity gradient: giant vs. 
faint and compact H n regions 

We now examine whether any selection effect can be responsi- 
ble for the difference between the steep gradient found in the 
early studies and the shallower gradient of this work. To this 
goal, we subdivided the sample of H n regions according to their 
projected size and surface brightness in the Ha emission-line. 
Then, we computed the intrinsic luminosity and the radius in an 
Ha emission-line calibrated map (courtesy of R. Walterbos) for 
each nebula of the whole sample. Defining the surface bright- 
ness (SB) as the ratio between the total flux and the area ex- 
pressed in arcsec 2 , we subdivided the sample according to their 
size and SB. Considering their size, we defined them as small if 
their radius R <15" (60 pc at the distance 840 kpc) and as large 
if R >15". Considering their surface brightness, we define them 
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Table 6. Giant H n regions with derived chemical abundances. 



Id 


Rgc 




R 


SB 




fknrl 
Ww 


n0~ 15 ero/rm 2 O 


(til cscc) 


(10 crg/cm s aicscc ) 


NGC595 


1.7 


13.0 


25 


16.1 


COOlAb 


1.9 


3.2 


20 


6.2 


MA2 


2.5 


4.5 


25 


5.6 


BCLMP691 


3.3 


2.4 


15 


8.2 


NGC604 


3.5 


36.6 


30 


31.4 


IC131 


3.9 


1.9 


15 


6.6 


BCLMP290 


4.2 


2.2 


15 


7.5 


NGC588 


5.5 


4.2 


20 


8.0 


IC132 


6.4 


3.5 


20 


6.8 



as bright if their surface brightness SB> 5.5xl0~ 19 erg cirT 2 s _I 
arcsec 1 , and as faint if SB is lower than this limit. The four 
combinations are allowed, i.e. H n regions can be small and ei- 
ther bright or faint, or large and again bright or faint. In Table|6] 
we show the galactocentric distance, Rcc> the Ha observed total 
flux, Ff{ a , the radius, R, and the SB, of the so-called giant re- 
gions, i.e. those with SB> 5.5xl0~ 19 erg cirT 2 s arcsec' 1 and 
R>15". 

To compare the total population with the giant H n regions, 
we show in Fig.|4]the oxygen abundances of the cumulative sam- 
ple, averaged in bins of 1 kpc each, together with the abundances 
of each single giant H n region. The O/H gradient of the H n re- 
gions in Table [6] computed with a weighted linear least-square 
fit, is 



12 + log(0/H) = -0.089(+0.023) R GC + 8.72(+0.09). 



(4) 



The gradient of the remaining sample is the same of given 
in Eq[TJ The giant regions show a significantly steeper gradi- 
ent, consistent with the gradients by Smith (119751 1. Kwitter & 
Aller (USD, Vflchez et al. ( fT988j l, and Garnett et al. JT9971 ). 
The question is whether this gradient is really significantly dif- 
ferent from the whole sample, and if this is the case, what is the 
reason of such different behaviour. 

Owing to the small number of giant H n regions (9 in to- 
tal), the uncertainty on the slope of their gradient is high. Thus 
it could still be in partial agreement, within the errors, with the 
larger sample, and their difference might stem from metallicity 
fluctuations in the ISM. On the other hand, the characteristics 
of the giant regions might be truly different from the average 
sample. For example large self -bound units are not destroyed by 
massive stars and thus retain their original structure and get con- 
tinuously enriched by SF. However, while this might be the case 
for the metallicity peak near the centre, it does predict that giant 
regions should have higher metallicity at all galactocentric radii, 
which is clearly not the case. Moreover, in a recent paper, Relano 
& Kennicutt (2009) studied the star formation in luminous H n 
regions in M33, which correspond mostly to our giant H n re- 
gions. They found that the observed UV and Ha luminosities are 
consistent with a young stellar population (3-4 Myr), born in an 
instantaneous burst. Thus the steeper gradient might result form 
a combination of a small statistics and of a metal self-enrichment 
effect in the giant region sample. 
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Fig. 4. The O/H radial gradient: the giant H n regions (empty 
triangles) and the complete sample averaged in bins, each 1 kpc 
wide (filled squares). The continuous line is the weighted mean 
least square fit of the giant H n region sample, while the dotted 
line refers to the complete sample. 



their radial gradient, but also their spatial distribution projected 
onto the disk. In Fig. [5] we show the two-dimensional metallicity 
distributions for M33 from H n regions and from PNe superim- 
posed to a contour map of the stellar mass distribution derived 
from the JHK image, a composition of the image of Regan & 
Vogel (119941) and the 2MASS image. The O/H abundances were 
averaged in bins of 0.8x0.8 kpc 2 . The white pixels indicate ar- 
eas where metallicity measurements are lacking: for H n regions 
they generally correspond to the interarm regions, while for PNe 
to spiral arms. 

The H n regions with the highest metallicity are not located 
in the optical centre of the galaxy (0,0 in the map), but rather lie 
at radius 1-2 kpc in the southeast direction. Also in the case of 
PNe, most of the metal rich PNe are located in the southern part 
of M33, from 2 to 4 kpc from the centre. However the lack of 
known PNe in the northern spiral arm at the same distance of the 
southern metal-rich PNe (because of the extended H n regions 
not allowing the identification of stellar emission-line sources) 
does not allow a complete 2D picture of their metal distribution 
around the central regions. 

To estimate the location of the off-centre metallicity peak 
for the H n region map, we divided its squared 10'xlO' region, 
centered at RA 1:33:50.9 dec 3 0:39:3 6 (M33 centre from the 
2MASS survey, Skrutskie et al. [20061 . with a 10x10 grid. We 
computed the radial O/H gradients varying the central position 
in the grid and then finding the one that minimizes the scatter 
of the gradient. We found an off-centre position at RA 1:33:59 
dec 30:33:35 (J2000.0), which corresponds to the location of the 
high-metallicity H n regions in the southern arm. The oxygen 
gradient measured from this central position of the whole H n 
region population, including also the central objects, is 



3.4. The 2-dimensional distribution of metals 

The usual way to study the metallicity distribution in disk galax- 
ies is to average it azimuthally, assuming that i) the centre of 
the galaxy coincides with the peak of the metallicity distribu- 
tion and if) the metallicity distribution is axially symmetric. The 
large number of metallicity measurements in M33, both from 
H n regions and from PNe, allowed us to reconstruct not only 



12 + log(0/H) = -0.021 (+0.007) R GC + 8.36(+0.03), (5) 

and it is shown in Fig. [6] The gradient of Eq.[5]is flatter and with 
a somewhat lower absolute dispersion than the one in EqQ] The 
reduction of the dispersion in the O/H gradient due to the dis- 
placement of the galaxy centre is not enough strong to confirm 
that the metallicity maximum corresponds to the real centre of 
the metallicity distribution. 
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Fig. 6. The O/H radial gradient computed with the centre located 
7 5 atRA 1:33:59 dec 30:33:35 (J2000.0) to minimize the dispersion 
in the slope of the radial gradient. 



tional centre in the galaxy; Hi) the selection criterion of H n re- 
gions for metallicity determination. 




1000 
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Fig. 5. The oxygen abundance maps of M33 (60'x 60'): H n re- 
gions (top) and PNe (bottom). The O/H abundances are aver- 
aged in bins of 0.8x0.8 kpc 2 . The colour-scale shows the oxy- 
gen abundance as indicated in the label. North at the top, east to 
the left. The optical centre of M33 is located at (0,0). The con- 
tour levels represent the stellar mass distribution derived from 
the JHK image ofM33. 



4. Why is the metallicity peak off-centre? 

In the following, we examine several possibilities to explain 
the presence of the off-centre metallicity maximum and the low 
metallicity measured in the central region: z) a local effect of 
ISM metallicity fluctuations; z'z) the lack of dominant gravita- 



4.1. Local metallicity fluctuations 

Simon & Rosolowsky (20083 have already noticed the non- 
axisymmetric distribution of H n region abundances and sug- 
gest that the material enriched by the most recent generation 
of star formation in the arm might not have been azimuthally 
mixed through the galaxy. A strong OB association located in 
the southern arm might be responsible for the enhancement en- 
richment at the location of the metallicity peak. Velocity shear 
is present in M33 even at small radii because of the slow rise 
in the rotation curve (Corbelli et. 2003). At the peak location, 
differential rotation will shear up the bubble of metals produced 
in the surrounding ISM in about 10 8 yrs. The timescale seems 
long enough to allow vigorous star formation at a particular lo- 
cation to enrich the ISM of metals well above the average value. 
However, the large dispersion in the metallicity around the peak 
location seems to rule out an inefficient azimuthal mixing or re- 
distribution of the metals. 

4.2. The lack of a gravitational centre 

The non-axysymmetric metallicity distribution might be re- 
lated to a general non-axisymmetric character of central re- 
gions of M33, noticed in the past by several authors. Colin & 
Athanassoula (1981) found that the young population displace- 
ment is located towards the southern side of M33 and amounts 
to approximately 2-3 ', i.e. 480-720 kpc. Using evidence of other 
asymmetries in the inner regions of M33, such as those present 
in the distribution of H i atomic gas, of H n regions, and in the 
kinematics, they proposed a bulge centre presently located in 
the northern part of the galaxy, which is rotating retrogradely 
around the barycenter of the galaxy. The analysis of infrared 
images (Minniti et al. [T993), however, seems to point out to a 
small bulge with a much smaller displacement the one advo- 
cated by Colin & Athanassoula (1981). A detailed analysis of 
the kinematics of the innermost regions of M33 by Corbelli & 
Walterbos (2007) confirms asymmetries in the stellar and gas ve- 
locities, which however seem more related to the presence of a 
weak bar. The exact galaxy centre is uncertain on scales of a few 
arcsec. Thus, even if M33 lacks of a dominant gravitational cen- 
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Table 7. H n regions in the central regions without direct 
r„measurement. 



ID 


RA 


Dec 


Rgc 


T e 


12 + log(0/H) 




J2000.0 


(kpc) 


(K) 


dex 




(1) 


(2) 


(3) 


(4) 


(5) 


BCLMP 93a 


1:33:52.6 


30:39:08 


0.23 


8700 


8.53" 


BCLMP 301h 


1:33:52.6 


30:45:03 


1.45 


9200 


8.60 


BCLMP 4a 


1:33:59.5 


30:35:47 


1.48 


9200 


8.33 


M33-SNR 64 


1:34:00.1 


30:42:19 


0.87 


8950 


8.40 



Vflchez et al. 1988 (CC93) adopted 7>6000 k obtaining 12 + log (O/H)=9.02 
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Fig. 7. The radial gradient of T e for the H n regions of the present 
sample and literature data. The continuous line is a mean least 
square fit. 

tre of M33 and the bright central cluster might migrate around 
it, it seems unlikely that the centre of the galaxy is off by several 
hundreds pc from where the bright cluster lie. The marginal gain 
in the dispersion in re-computing the metallicity gradient from 
an off-centre position (see previous section) confirms that this 
hypothesis seems unlikely. 

4.3. Selection criterion 

That the average metallicity at the centre seems lower than at 
1 .5 kpc is hard to explain in the framework of an inside-out disk 
formation scenario. We now discuss the possibility than in the 
central regions the metallicity might be higher than reported in 
this paper because of a bias in the H n region selection. As ex- 
plained in Sect. 3.1 the inclusion of H n regions in our sample 
requires determining of the electron temperature through the de- 
tection of the faint oxygen auroral line. As the metallicity in- 
creases the line becomes so faint as to be detectable only in 
bright complexes. The centre of M33 lacks of vigorous star- 
forming sites, so the most cooler, metal rich H n regions have the 
oxygen auroral line below the detection threshold. To prove that 
this might be the case we searched the literature for the existing 
H ii region spectra inside 1 .5 kpc radius, which were not included 
in our database because of the undetectable [O m] 4363 line. We 
found 4 H ii regions in the database of Magrini et al. (|2007a) 
for which optical spectroscopy is available but no detection of 
temperature diagnostic lines. Their names, coordinates, galac- 
tocentric distances, assumed electron temperatures, and oxygen 
abundances from M07a fluxes are shown in TableQ 

In Fig|7]we plot the relationship between the electron tem- 
perature and the galactocentric distance for the complete H n 
region sample. The weighted mean least square fit gives a rela- 
tionship between the two quantities 

T e = (410 ± 80) x R cc + 8600 + 320 (6) 

where T e is expressed in K and Rcc in kpc. We thus derived 
T e for the H ii regions of Table [7] adopting it in their chemical 
abundance calculations. 

Using the intensities by M07a, complemented with upper 
limit measurements of the [O n] 7320-7330 A when not avail- 
able in the original paper, we roughly estimated the oxygen 
abundance of the four H n regions located within 1.5 kpc. Their 
location in the radial gradient is shown in Fig. [8] The sources 
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Fig. 8. The O/H radial gradient: the symbols for the H n regions 
with a positive detection of T e diagnostic lines are as in Fig. 1, 
while the filled squares are the regions of Table|7]for which the 
T e is extrapolated from Eq. [6] The continuous line is the fit of 
Eq. 4. 



clearly lie above the average metallicity determined in the centre 
of our database. It is therefore likely that the adopted criterion 
for source selection, based on the positive detection of lines for 
determining the electron temperature, might be responsible for 
the low metallicity derived in the centre. 

In the 1.5 kpc central region, there are 8 H n regions 
with measured electron temperature having oxygen abundances 
from ~8 to -8.5 (see Fig. 1), namely B0043b (O/H 8.214), 
B0029 (8.211), B0038b (8.391), B0016 (7.98), B0079c (8.386), 
B0027b (8.37), B0033b (8.399), B0090 (8.506). Their metal 
abundance and dispersion are consistent with the hypothesis that 
we are missing several high-metallicity regions within 1.5 kpc 
of the centre since they only seem to trace the low metallicity 
side of the distribution present at each given radius. Similarly, 
this bias partially explains the different gradient of the giant H n 
regions: we only include metal-rich H n regions (with metallic- 
ity above a critical value) in the sample if they are very luminous 
because only these show detectable temperature diagnostic lines. 

5. The chemical evolution of M33 

The metallicity gradient derived in Sect B.ll from our sample of 
H ii regions characterizes the ISM composition in M33 at the 
present time. Together with the metallicity gradient from the PN 
population (M09), these data allow setting new constraints on 
current models of galactic chemical evolution. 

The model of M07b, specifically designed for M33, is able 
to predict the radial distribution of molecular gas, atomic gas, 
stars, SFR, and the time evolution of the metallicity gradient. In 
the following, we discuss the modifications needed to reproduce 
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the metallicity gradient of H n regions and PNe derived in this 
work and in M09. 



5.1. A revised model of chemical evolution 

The multiphase model adopted by M07b follows the formation 
and destruction of diffuse gas, clouds, and stars, by means of the 
simple parameterizations of physical processes (e.g. Ferrini et 
al. I1992I I. In particular, the SFR is the result of two processes: 
cloud-cloud interactions (the dominant process) and star forma- 
tion induced by the interaction of massive stars on molecular 
clouds. At variance with other models, the relationship between 
the SFR surface density, Zsfr, and the molecular gas surface 
density, E# 2 , (the so-called Schmidt law) is thus a by-product 
of the model and is not assumed. In general, the relation be- 
tween the surface density of total gas and SFR has a slope of 
1.4+0.1, but the slope can vary from galaxy to galaxy (Kennicutt 
1998a). In the particular case of M33, a tight correlation exists 
between the SFR, measured from the FUV emission, and the sur- 
face density of molecular gas has a well-defined slope (Verley et 
al. 2009), 



Stars form by cloud-cloud interactions at a rate H and by the 
interactions of massive stars with clouds at a rate a, 



= AS 



1.1±0.1 



(7) 



According to M07b, the best model for M33 (the so-called 
accretion model, with almost constant infall) suggests a long 
lasting phase of disk formation resulting from a continuous ac- 
cretion of intergalactic medium during the galaxy lifetime. We 
refer to M07b for the general description of the model and of 
the adopted parameters for M33. Here we concisely describe the 
model and only the updated processes and equations in detail. 

The galaxy is divided into N coaxial cylindrical annuli with 
inner and outer galactocentric radii /?, (z = l,N) and Rm, 
respectively, mean radius /?;+i/2 = (Ri + Ri+\)/2 and height 
/i(R, + i/2). Each annulus is divided into two zones, the halo and 
the disk, each made of diffuse gas g, clouds c, stars s, and stellar 
remnants r. The halo component includes the primordial bary- 
onic halo but also the material accreted from the intergalactic 
medium. 

At time t — 0, all the baryonic mass of the galaxy is in the 
form of diffuse gas in the halo. At later times, the mass frac- 
tion in the various components is modified by several conversion 
processes: diffuse gas of the halo falls into the disk, diffuse gas 
is converted into clouds, clouds collapse to form stars and are 
disrupted by massive stars, and stars evolve into remnants and 
return a fraction of their mass to the diffuse gas. In this frame- 
work, each annulus evolves independently (i.e. without radial 
mass flows) keeping its total mass (halo + disk) fixed from t — 
tof ga i = 13.6 Gyr. 

The disk of mass Mo(t) is formed by continuous infall from 
the halo of mass M#(f) at a rate 



dM L 
dt 



= fM H , 



(8) 



where / is a coefficient proportional to the inverse of the infall 
timescale. Clouds condense out of diffuse gas at a rate jj and are 
disrupted by cloud-cloud collisions at a rate H' , 



dt 



mV 2 - h'mI, 



(9) 



where M g (t) is the mass fraction (with respect to the total bary- 
onic mass of the galaxy) of diffuse gas, and M c (t) was defined 
above. 



dM s 



HM 2 C + aM s Mc ~ DM S , 



(10) 



where M s (t) is the mass fraction in stars, D the stellar death rate 
and M c (t) as above. All rate coefficients of the model are as- 
sumed to be independent of time but are functions of the galac- 
tocentric radius Rqc (cf. M07b for their values and radial behav- 
iors). 

5.2. The choice of the stellar yields and the IMF 

The model results are sensitive to the assumed stellar yields. For 
low- and intermediate-mass stars (M < 8 M Q ), we used the 
yields by Gavilan et al. d2005l . The yields of Marigo d200T1> 
give comparable results without any appreciable difference in 
the computed gradients of chemical elements produced by in- 
termediate mass stars, such as N with respect to the gradients 
computed with the yields by Gavilan et al. (2005). For stars in 
the mass range 8 M Q < M < 35 M we adopt the yields by 
Chieffi & Limongi (2004). The yields of more massive stars are 
affected by the considerable uncertainties associated with differ- 
ent assumptions about the modelling of processes like convec- 
tion, semi-convection, overshooting, and mass loss. Other dif- 
ficulties arise from the simulation of the supernova explosion 
and the possible fallback after the explosion, that strongly in- 
fluences the production of iron-peak elements. It is not surpris- 
ing then that the results of different authors (e.g. Arnett 1995; 
W oosley & Weaver 119951 Thielemann et al. 119961 Aubert et 
al. [1996) disagree in some cases by orders of magnitude. In 
our models, we estimate the yields of stars in the mass range 
35 M e < M < 100 M e by linear extrapolation of the yields in 
the mass range 8 M < M < 35 M Q . 

Another important ingredient in the chemical evolution 
model is the initial mass function (IMF). Several works support 
the idea that the IMF is universal in space and constant in time 
(Wvse [T997l Scalo [T9981 Kroupa 120021 . apart from local fluc- 
tuations. There are several parameterizations of the IMF (see 
e.g. Romano et al. 2005 for a complete review), starting from, 
e.g. , Salp eter ( fl955l ), Ti nsley ( [19801 ), Scalo | fl5§6l , Kroupaet 
al. JT9931 Ferrini et al. (fl99l)l l. Scalo ( [19981 . Chabrier ( T20031) . 
In the following, we test the possibility that the observed flat gra- 
dients can be explained in terms of a non-standard IMF. In fact, 
the magnitude and the slope of chemical abundance gradients 
are related to the number of stars in each mass range, and so to 
the IMF. The goal is to reproduce, if possible, the flat gradient 
supported by recent observations by only varying the IMF. 

In Fig. [9] we compare the oxygen and nitrogen gradients 
of H 11 regions with present-day abundance profiles from the 
M07b model. The choice of the IMF does not affect the slope 
of O/H and N/H gradients, but simply shifts the abundance pro- 
files to higher (or lower) values according to the amount of 
stars that produce oxygen (massive stars) or nitrogen (low- and 
intermediate-mass stars). With the adopted stellar yields, the best 
fits of the metallicity distribution of M33 are obtained with the 
IMF by Ferrini et al. (1 19901 1 and Scalo dl9861 l. For the revised 
model we thus adopt the parameterization of Ferrini et al. (1990). 
The different slopes of the O/H and N/H gradients are repro- 
duced fairly well by the chemical evolution model as a natural 
consequence of the different mass ranges of stellar production, 
hence timescales, of these two chemical elements (see Sect l3.2l 
for a comparison with the measured gradients). 
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Fig. 9. The oxygen (top panel) and nitrogen (bottom panel) ra- 
dial gradients of M07b with several parameterization of the IMF: 
Salpeter (1955) (continuous green line), Tinsley (1 19801 1 (dotted 
red line), Scalo (1 19861 1 (dashed magenta line), Scalo 1998) (long 
dashed yellow line), Ferrini et al. ( 1990) (long dash-dotted blue 
line), Chabrier (2003) (dash-dotted black line). 



5.3. Comparison with the observations 

The observational constraints to the model are those described in 
M07b, complemented with the O/H and N/H radial gradients of 
H n regions and PNe from this work and M09, respectively, and 
the radial profile of the SFR determined by Verley et al. (2009) 
from far ultraviolet (FUV) observations corrected for extinction. 
For H ii regions we used the gradient derived from the whole 
population (see Eq. 4), without any distinction in terms of size 
and brightness, but excluding the central 1 kpc. Giant H n re- 
gions might not be representative of the current ISM abundance 
due owing their possible chemical self-enrichment. 

In Figs. [10] and [TJJ we show the radial surface den- 
sity of molecular gas (from the single dish observations of 
Corbelli 2003 averaged over 1 kpc bins) and the SFR as a func- 
tion of the surface density of molecular gas (Verley et al. 2009), 
respectively. Both figures show the predictions of the model of 
M07b at t = fgai. Clearly, the model of M07b in its original 
formulation is unable to reproduce these constraints. We have 
therefore considered other parameterizations of the star forma- 
tion process by varying the radial dependence of the SFE, repre- 
sented by the coefficient H. 

Our experiments show that, to reproduce the observed be- 
haviour of the radial gas distribution and the observed Schmidt 
law, is necessary to increase the efficiency of star formation at 
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Fig. 10. The SFR derived from the UV emission corrected for 
extinction by Verley et al. (2009) vs. the surface density Schmidt 
law: filled circles (blue) are the molecular gas averaged in bins 
1 kpc each (Corbelli 120031 curves represent the model at 0.5 
Gyr from the disk formation (dotted curve), 2 (dashed curve), 3 
(long-dashed curve), 5 (dot-dashed curve), 8 (long dash-dotted 
curve), 12 (long-short dashed curve), and at 13.6 Gyr (solid red 
curve); the solid green curve is the model by M07b at 13.6 Gyr. 
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Fig. 11. The radial surface density of molecular gas: filled circles 
(blue) are the observed molecular gas averaged in bins 1 kpc 
wide. Model curves have the same symbols as in Fig.fTOl 

large radii. This can be accomplished in many ways. In the pre- 
vious version of the model (M07b) H decreased with Rqc to con- 
sider a geometrical correction resulting from how large galacto- 
centric distances correspond to larger volumes (see Ferrini et al. 
II 9941 ) . In the present work we assume that H is constant with 
radius, thus implying that the star formation efficiency increases 
linearly with galactocentric radius. In Sect. I5.4l we describe the 
observational evidence in support this assumption in M33. 

Assuming that H is spatially constant, we obtain the results 
shown in Figs. [U and [TT] at t = 0.5, 2, 3, 5, 8, 12, 13.6 Gyr. 
The relationship between molecular gas and SFR predicted by 
the modified M07b model is more or less constant with time and 
in good agreement with the data. The predicted Schmidt law (at 
the present time) has an average exponent ~ 1 .2, in agreement 
with the observations, whereas the original M07b model pro- 
duces a Schmidt law with a higher exponent, ~ 2.2. The better 
agreement with the observations obtained with the revised M07b 
model suggests that M33 is more efficient in forming stars than 
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normal local Universe spiral galaxies, in particular in its outer 
regions. 



5.4. Flat gradients and SF efficiency 

Is there evidence of flat gradients in other galaxies? How can 
they be explained? A flat gradient in the outer regions of our 
Galaxy has been observed by several authors (e.g., Yong et 
al. [20031 Carraro et al. [20071 Sestito et al. 120071 using differ- 
ent metallicity tracers (e.g., Cepheids, open clusters). However, 
other tracers, such as PNe, do not show this behaviour, indicating 
flat gradients across the whole disk (cf. Stanghellini et al. 2006 
Perinotto & Morbidelli 2006 ). The outer Galactic plateau might 
be a phenomenon similar to the flat gradient of M33. Since M33 
is less massive than the MW, its halo collapse phase, responsible 
for the steep gradient in the inner regions (Rcc < 1 1-12 kpc) of 
our Galaxy (cf. Magrini et al. I2009bl l. is less marked. Thus in 
M33 the difference between the inner and outer gradients is less 
evident than in the MW. The flat metallicity distribution of the 
MW at large radii has been explained with several chemical evo- 
lution models, among them those by Andrievsky et al. (120041) . 
Chiappini et al. (120011) model C, and Magrini et al. d2009bl ). 
Andrievsky et al. ([2004) explain the flat metallicity distribution 
beyond 11-12 kpc when assuming that the SFR is a combination 
of two components: one proportional to the gas surface density, 
and the other depending on the relative velocity of the interstel- 
lar gas and spiral arms. With these assumptions they explain why 
the breaking point in the slope of the gradient and the consequent 
outer flattening, occur around the co-rotational radius. Chiappini 
et al. ( 120011) assume two main accretion episodes in the lifetime 
of the Galaxy, the first one forming the halo and bulge and the 
second one forming the thin disk. Their model C is the one able 
to reproduce a flat gradient in the external regions, assuming that 
there is no threshold in the gas density during the halo/thick- 
disk phase, and thus allowing the formation of the outer plateau 
from infalling gas enriched in the halo. Finally, Magrini et al. 
(2009b) reproduce the Galactic gradient thanks to the radial de- 
pendence of the infall rate (exponentially decreasing with radius) 
and with the radial dependence of the star and cloud formation 
processes. To reproduce a completely flat gradient in the outer 
regions would, however, require additional accretion of gas uni- 
formly falling onto the disk, which would result in inconsistent 
behaviour of the current SFR. 

Why do we use a radially increasing SFE to reproduce the 
flat metallicity gradient of M33? M33 in general is quite differ- 
ent from large spiral galaxies in terms of SF. A comparison of 
the SFR to the H2 mass shows that M33, like the intermediate 
redshift galaxies, has a significantly higher SFE than large lo- 
cal Universe spirals. There is also observational evidence that 
the S FE var ies with radius in M33. Kennicutt ( 1998b), Wong & 
Blitz (2002 ), and Murgia et al. d2002) find a molecular gas de- 
pletion timescale (inversely proportional to the SFE) that varies 
radially, decreasing by a factor ~2-3. Also Gardan et al. (2007) 
from the CO to SFR ratio also find a dependence of the depletion 
timescale on radius, decreasing of a factor ~2 over 4 kpc. These 
observations prompted us to model the chemical evolution of 
M33 taking SF processes depending on radius into account. The 
general assumption of our model is that the SF is driven by cloud 
collisions. Increasing the efficiency of this process with radius 
does not imply that cloud-cloud collisions are more efficient in 
the outer regions, but that additional processes may contribute to 
the SF in the peripheral regions. 



5.5. Implications for the metallicity gradient and its evolution 

The introduction of more realistic SF process in M33 has signif- 
icant consequences on its metallicity gradient and evolution. As 
described by M09, oxygen modification (destruction or creation) 
does not occur in the PN population of M33. Thus oxygen, the 
best-determined element in nebular optical spectroscopy, can be 
safely used as a tracer of the ISM composition at the epoch of 
the formation of PN progenitors. Most PNe in M33 are older 
than 0.3 Gyr, and probably much older, up to an age of 10 Gyr. 

The metallicity gradient of PNe, including only those with 
progenitor stars older than 0.3 Gyr and excluding PNe located in 
the first kpc from the centre, is 

12 + log(0/H) = -0.040(+0.014) R GC + 8.43(+0.06). (11) 

As discussed above, the whole sample of H 11 regions is repre- 
sentative of the current ISM composition. The H 11 region and 
PN gradients in the same radial region, 1 kpc < Rcc 55 8 kpc, are 
indistinguishable within the uncertainties in their slopes, with a 
slightly translation, ~0. 1 dex, towards higher metallicity in the 
H n region sample. This means that very little evolution of metal- 
licity has occurred in the past few Gyr. 

In Figure [12] we show the observed metallicity gradient of 
M33 together with the results of the present model (top panel) 
and of the original model of M07b (bottom panel). For a better 
comparison, the oxygen abundances of the samples of H n re- 
gions and disk PNe have been averaged over bins of 1 kpc each. 
This allows us to highlight the translation towards the lower 
metallicity of the PNe. The curves correspond to the present and 
to 5 Gyr ago (8.6 Gyr from the disk formation, assuming an age 
of 13.6 Gyr for M33), i.e., approximately the average period of 
formation of the PN progenitors, the present model predicts a 
higher SF in the outer regions, reproducing the observed flatness 
of the oxygen gradient at large radii better than the M07b model. 
In addition, the evolution of the metallicity gradient predicted by 
the new model is consistent with the observations: a negligible 
change in the slope and a small translation to higher metallic - 
ities. Finally, we see that the assumed dependence of the SFE 
with galactocentric radius, required to reproduce the metallic- 
ity evolution of M33, takes place always as in M07b with a slow 
building up of the disk of M33 by accretion from an intergalactic 
medium or halo gas. 

For completeness we also show the N/H radial gradient of 
the sample of H 11 regions. For this element the temporal evolu- 
tion cannot be inferred using PN abundances, since PNe modify 
their nitrogen composition during their lifetime. The nitrogen 
abundances have been averaged over bins of 1 kpc each. The 
observed N/H gradient is significantly steeper than the oxygen 
gradient, as predicted by the model. 

The time evolution of the metallicity gradient is strictly cor- 
related to the inside-out growth of the disk. In a recent work, 
Williams et al. (2009) presented the resolved photometry of four 
fields located at different galactocentric distances, from 0.9 to 
6.1 kpc. Their photometry provides a detailed census of stellar 
populations and their ages at different radii. They find that the 
percentage of the stellar mass formed prior to z- 1 (~8 Gyr ago) 
changes from ~71% in the innermost field at 0.9 kpc from the 
centre to ~16% in the outermost SF region at 6.1 kpc. We thus 
compare these results with the cumulative SF predicted by our 
model at different ages. The agreement is good in the central 
field, where model and observations both have about 70% of the 
stellar mass already formed 8 Gyr ago. Also at 2.5 kpc, corre- 
sponding to their second field, the observed stellar mass formed 
prior to z- 1 is about 50% and the model results give 47%. In 
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Fig. 12. The O/H radial gradient evolution. Top panel: filled 
circles (blue) are the H n region oxygen abundances from the 
dataset described in Sectf2] averaged over radial bins, each 1 kpc 
wide; filled squares (magenta) are averaged non-Type I PN oxy- 
gen abundances by M09. The gradient 5 Gyr ago (dashed curve) 
and now (continuous curve) as predicted by the present model. 
Bottom panel: the same as in the top panel but for the accretion 
model by M07b. 

the outermost fields, at 4.3 and 6.1 kpc, the agreement is not as 
good, since the observations show that 20% and 16% of stars 
formed before 8 Gyr ago, while the model predicts ~44% for 
both fields. With the model by M07b, the situation is almost un- 
varied in the inner regions, while in the outer regions there is 
a slightly lower percentage of stars formed before 8 Gyr ago, 
~33%. The recent evolution of the outermost SF regions of M33 
can be explained by an accretion of material in the peripheral re- 
gions at recent epochs. In our model, the scalelength of the disk 
is constant, and the accretion in the outer regions is a continuous 
process. Thus, we are able to reproduce the integrated SF at any 
radius and its by-product, the metallicity, but not to reproduce its 
temporal behaviour if it is dominated by stochastic events, such 
as, e.g., sporadic massive accretion of gas and tidal interactions. 
Finally, outside the SF area of the disk, the age of populations of 
the outer disk/halo show an age increase with radius (Barker et 
al. 2007a, 2007b). This agrees well with our model results: there 
is a slightly larger number of old stars beyond 8 kpc than in the 
peripheral part of the SF disk. 

6. Summary and conclusions 

We have studied the chemical evolution of M33 by means of new 
spectroscopic observations of H n regions, together with litera- 
ture data of H n regions and PNe. We derived the radial oxygen 
gradient in the range from 1 to ~8 kpc. Its slope is -0.044± 0.009 
dex kpc 1 . We excluded the central 1 kpc region from the gra- 
dient because of its low metallicity. In fact, the 2D metallicity 
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Fig. 13. The N/H radial gradient. Top panel: filled circles (blue) 
are the H n region nitrogen abundances from the dataset de- 
scribed in Sect|2] averaged over radial bins, each 1 kpc wide. The 
gradient at the present time as predicted by the present model is 
the continuous curve. Bottom panel: the same as in the top panel 
but for the accretion model by M07b. 

map is off-centre, with a peak in the southern arm at 1-2 kpc 
from the centre. We measured the highest metallicity gradient in 
bright regions which are not located in the centre of M33. We 
explained this effect with a bias in the measurement of nebular 
abundances towards low metallicity in the central regions and/or 
with SF bursts which had not had time to mix azimuthally. At 
a galactocentric distance of about 1.5 kpc, the spread of metal 
abundances is much greater than the dispersion along the aver- 
age gradient, and it is possibly related to the presence of a bar. 

We analysed and discussed the metallicity gradient of giant 
H ii regions, i.e. the largest regions with high surface brightness. 
We find it steeper than the average of the sample gradient, in 
agreement with the early spectroscopic studies of H n regions in 
M33, based mainly on the brightest objects. 

We compared the metallicity gradients of H n regions and 
of PNe, obtained with the same set of observations and analysis 
techniques. We find a substantially unchanged slope of the gra- 
dient, and an overall increase in metallicity with time. We can 
explain the slow evolution of the metallicity gradient from the 
present to the birth of the PN progenitors with a chemical evolu- 
tion accretion model, as in M07b, if the SFE is higher at larger 
radii. 

Appendix A: Comparison with the literature 

To assess the quality of the MMT observations we compared the 
oxygen abundance in the control sample with those available in 
the literature. Oxygen is in fact the most widely used element 
in emission-line optical spectroscopy because it gives accurate 
results for determining the metallicity. This is because /) many 
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ionization stages are present in the optical range ([O i], [O n], 
[O in]), and the corrections for unseen ionization stages are of- 
ten not needed ; ii) a direct measurement of its electron temper- 
atures, T e [0 III] and T e [0 II], is possible. For each source in the 
control sample, we show in TablelAlthe oxygen abundance deter- 
mined by different authors (Rosolowsky & Simon 2008 Magrini 
et al. l2007al Crockett et al. l2006l Vilchez et al. [19881 Kwitter & 
Aller 1981 ) and in this paper. In columns 2 to 6 we report the lit- 
erature O/H and our measurements in column 7. For H n regions 
with multiple determinations we computed the weighted average 
of all measurements. The agreement is reasonably good for all 
sources, since our values are consistent, within the errors, with 
those available in the literature. The average O/H abundance of 
the 14 H n regions agrees. The average dispersion between our 
values and previous determinations is only 0.09 dex. 
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